rm(list=ls())

#Ns=c(30,50,100,200,300,400,500,600,800,1000)
#Ks=c(1,2,3,4,5,6,7,8,9,10)

# these Ns and K's are used in SM's fig2.
Ns=c(100,200,300,400,500)
Ks=c(1,2,3,5,8,10) # corresponding to k=2,3,4,6,9,11 in SM (they count T)


datsim=read.table("NK.G.out")

dim(datsim)

Dat={}

for(n in 1:length(Ns)){
  for(k in 1:length(Ks)){
    dat=datsim[datsim[,1]==Ns[n]&datsim[,2]==Ks[k],]
    Dat=rbind(Dat,apply(dat,2,mean))
  }
}

dim(Dat)

N=Dat[,1]
K=Dat[,2]
nearest=Dat[,4]  # mean(dist from cf to nearest match);

nearestk<-matrix(0,length(Ns), length(Ks))
for(k in 1:length(Ks)){
nearestk[,k]=nearest[K==Ks[k]]
}

pdf(file="NK.G.pdf")

op <- par(mfrow = c(1, 1), cex=1.2)

matplot(Ns,nearestk,type="l",xlab="Number of Data Points",ylab="Average Distance to the Nearest Match")
text(150,.06,expression(k==2),cex=.8)
text(150,.16,expression(k==3),cex=.8)
text(150,.22,expression(k==4),cex=.8)
text(150,.28,expression(k==6),cex=.8)
text(150,.33,expression(k==9),cex=.8) 
text(150,.38,expression(k==11),cex=.8)

par(op)

dev.off()
q()
